%function mesh(N,X,Y,Nx,Ny,r)
h=1/8;
Lx=-1;
By=-1;
Rx=1;
Ty=1;
r=0.5;
%%%%%%%%%%%%%%%%%%%%%%%
Nx=(Rx-Lx)/h+1;%the number of node in  x direction
Ny=(Ty-By)/h+1;%the number of node in  y direction
 N=Nx*Ny;     %the number of grid point
 A=[0,0];  % coordinate of the circle  center
%%%%%%%%%%the node coordinate in Cartesian mesh%%%%%%%%
    X=Lx:h:Rx;
    Y=By:h:Ty;
 node=zeros(N,2);%the number of grid points
for  i=1:Ny
   coordinates((i-1)*Nx+1:i*Nx,1)=X(1:Nx);
   coordinates((i-1)*Nx+1:i*Nx,2)=Y(i);
end
% 
% %%%%%%triangle elememt %%%%%%%
Ntri=2*(Nx-1)*(Ny-1);%the number of triangle
elements3=zeros(Ntri,3);%node script of every triangle
for i=1:Ny-1
    elements3((i-1)*(Nx-1)+1:i*(Nx-1),1)=(i-1)*Nx+2:i*Nx;
    elements3((i-1)*(Nx-1)+1:i*(Nx-1),2)=i*Nx+2:(i+1)*Nx;
    elements3((i-1)*(Nx-1)+1:i*(Nx-1),3)=(i-1)*Nx+1:i*Nx-1;
    elements3((i+Nx-2)*(Nx-1)+1:(i+Nx-1)*(Nx-1),1)=i*Nx+1:(i+1)*Nx-1; 
    elements3((i+Nx-2)*(Nx-1)+1:(i+Nx-1)*(Nx-1),2)=(i-1)*Nx+1:i*Nx-1;
    elements3((i+Nx-2)*(Nx-1)+1:(i+Nx-1)*(Nx-1),3)=i*Nx+2:(i+1)*Nx; 
end

%%%%%%%%%the boundary%%%%%%%%%
Dirichlet=zeros(2*(Nx+Ny-2),2);
%%%%%%%left%%%%%%%%
Dirichlet(1:Ny-1,1)=1:Nx:(Ny-2)*Nx+1;
Dirichlet(1:Ny-1,2)=Nx+1:Nx:(Ny-1)*Nx+1;
%%%%%%%%%%%top%%%%%%%%%%
Dirichlet(Ny:Ny+Nx-2,1)=(Ny-1)*Nx+1:Nx*Ny-1;
Dirichlet(Ny:Ny+Nx-2,2)=(Ny-1)*Nx+2:Nx*Ny;
%%%%%%%%%%right%%%%%%%%%%%
Dirichlet(Ny+Nx-1:2*Ny+Nx-3,1)=Nx*Ny:-Nx:2*Nx;
Dirichlet(Ny+Nx-1:2*Ny+Nx-3,2)=Nx*(Ny-1):-Nx:Nx;
%%%%%%%%%%%bottom%%%%%%%%%%
Dirichlet(2*Ny+Nx-2:2*(Nx+Ny-2),1)=Nx:-1:2;
Dirichlet(2*Ny+Nx-2:2*(Nx+Ny-2),2)=Nx-1:-1:1;

% elem = label(node,elem);
trisurf(elements3,coordinates(:,1),coordinates(:,2),ones(size(coordinates,1),1))
view(2), axis equal, axis on,    drawnow;